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Abstract 

We argue that large fluctuations in the saturation momentum are necessary to 
explain the ATLAS and ALICE data on pA collisions measured at the LHC. 
Using a form for the distribution of fluctuations motivated by theoretical 
studies of the non-linear evolution equations for the Color Glass Condensate, 
we And a remarkably good agreement between theory and the measured 
distributions. If the saturation momentum fluctuates, we argue that the 
cross section for a proton probe should also fluctuate, consistent with previous 
observations. 


1. Introduction 

In a recent paper by Bzdak and Skokov [1] , it was argued that the num¬ 
ber of participants dependence of the multiplicity as a function of centrality 
provided a crucial test which might discriminate between the description of 
heavy ion collisions provided by the Color Glass Condensate [2] and that of 
the wounded nucleon modeljSj H]. Their basic idea is very simple: In the 
simplest version of the wounded nucleon model, the rapidity dependence of 
the multiplicity distribution for a collision of a projectile nucleon against 
^part ~ 1 nucleons in the target is a simple triangular distribution going from 
a number proportional the number or participants for rapidities close to that 
of the nucleus to a number of order 1 at rapidities corresponding to the beam 
nucleon. This is shown in Fig. [T] 
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Figure 1: An illustration of the multiplicity distribution of produced particles in pA col¬ 
lisions as a function of rapidity in the simplest version of the wounded nucleon model. 

On the other hand, the simplest version of the Color Glass Condensate 
model is that the rapidity distribution in the fragmentation of the target 
nucleus is of order iVpart but elsewhere of order 1 up to logarithms of the 
number of participants. This is because for a probe which has a Lorentz 
gamma factor of order 7 > relative to the nucleus, the nucleons in 
the nucleus act coherently, and for the typical transverse momentum scale 
associated with particle production, the nucleus appears as a black disk. The 
multiplicity is controlled by the number of partons in the proton evolved 
to the saturation scale of the nucleus. In fact there is a mild logarithmic 
correction to this, so the multiplicity scales as In iVpart [3 E] . Schematically 
the Color Glass Condensate description is shown in Fig. 

dN ~ Npart 



Figure 2: An illustration of the multiplicity distribution of produced particles in pA col¬ 
lisions as a function of rapidity in the simplest version of a model for the Color Glass 
Condensate. 

The number of participants dependence of the pseudorapidity distribu¬ 
tions in pA collisions at the LHC energy has been recently measured by the 
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ATLAS [7] and ALICE [S] collaborations. In ATLAS, they extract the num¬ 
ber of participant dependence of the integrated multiplicity under a variety of 
different assumptions concerning the cross-section for the proton penetrating 
through the nucleus. These correspond to a hxed Glauber cross section for 
the proton and to two different parameterizations of the Glauber-Gribov de¬ 
scription PEU]. The Glauber-Gribov model allows the proton cross-section 
to vary by Gausssian fluctuations. Similar procedures have been also applied 
by the ALIGE experiment. 

In the analysis of the next section, we compare expectations from the 
naive versions of the wounded nucleon model and that of the Golor Glass 
Gondensate to the experimental data. This analysis consists of comparing 
the Apart dependence predicted by the two models as a function of rapidity to 
that of the experimental data. Neither model gives an acceptable description 
of the data for any extraction of the number of participant dependence. 

How might one resolve this impasse? A hint is suggested by the extraction 
of the relationship between number of participants and multiplicity used to 
determine centrality. One might argue that in the fragmentation region of 
the nucleus, all models agree that the multiplicity dependence is linear in 
the number of participants. If so, the closest we come to the fragmentation 
region with the data is for 2 < r] < 2.7. While this is not forward enough 
to really be in the nuclear fragmentation region, we see that the data favors 
large fluctuations in the proton cross section. This suggests there might be 
fluctuations in the proton saturation momentum. 

Before we discuss fluctuations of the saturation scale, let us remark that 
there exists a number of theoretical papers that took an effort to describe 
the pA data within the framework of the different versions of the GGG-based 
models. A complete list of such models can be found in a review paper 
m that was, however, published before the hrst LHG run. More recently 
successful description of the pA data has been presented in Refs. 

[T5] and references therein. 

One expects a relationship between the cross section and the saturation 
momentum. If we require that the local density of gluons, denoted as A^, in 
an impact parameter model of the proton factorizes as proportional to the 
number of gluons per unit area on average times an impact parameter prohle, 
which for simplicity can be assumed to be exponential, then 

A2 = (1) 
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Requiring that the cross-section corresponds to the radius at which this den¬ 
sity has some hxed value Ag, perhaps we conclude that the radius of 

the proton grows as 

~ i ln(Q.VAa (2) 

jj, 

Therefore the cross section a oc will fluctuate if the saturation 

momentum does. There are of course many weaknesses of this analysis, 
perhaps the most important one is that we do not really know the impact 
parameter distribution for matter inside the proton, so that the exponential 
form, while valid at very large radii, may be misleading. 

The emergence of the saturation scale and geometrical scaling in the the¬ 
ory of the Color Glass Condensate has been introduced in Refs. [IS1II7|. It 
is based on the fact that the nonlinear evolution equations in small Bjorken 
X regime [iHl HH] possess traveling wave solutions [MllIT]. It has been ar¬ 
gued, however, in Ref. [23] and expanded in Ref. |2l| that the scaling CGC 
solutions have exponential tails which are unphysical, and therefore have to 
be modihed. This modihcation is based on the analogy with the statistical 
physics, and has clear interpretation in the dipole model, where the small 
tail of the amplitude is essentially counting a number of interacting dipoles. 
This discrete nature of the tail introduces stochastic diffusion that affects the 
shape of the entire amplitude, which - event by event - satishes geometrical 
scaling, with the saturation scale that is, however, somewhat different than 
the one obtained in the mean held approximation. It has been argued that 
these huctuations are Gaussian in the logarithm of the saturation momentum 
j23j . and controlled by a probability distribution: 

1 ( (\nQllQl-\nQllQlf 

= -- 

where 

P = (4) 

Here is the proton saturation momentum huctuating around its loga¬ 
rithmic average denoted as liiQ^ (with Qg ~ 1 GeV^/c^ being an arbitrary 
momentum scale, which cancels out in (|^). 
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It has been furthermore argued in |22] that the width a of the distribution 
(|^ is growing with the energy of the collision and depends upon the rapidity 
where measurements are made. One has to stress, however, that this result 
is only asymptotic, and - as the authors of Ref. |23] admit explicitly - there 
is no general analytical proof of this result and its status is still that of a 
conjecture. Phenomenologically scaling violations corresponding ro 

a ^ y/y with y = log(l/a;). (5) 

have been searched for in the DIS data with negative result [201 ED E2] • 
There are two possible ways out from this situation. Either ([^ holds down 
to the present energies (or presently accessible Bjorken x’s), in which case 
the distribution (j^ collapses to 6{\nQl/Qp), or result ^ is only asymptotic, 
and for low energies a tends to a constant. As we shall ^ow in the following, 
in such a case geometrical scaling would still hold. 

To this end let us consider a DIS amplitude T that exhibits geometrical 
scaling with respect to the fluctuating saturation momentum Qs{x): 

TiQ^Ql) = TipQ - p) (6) 

where pg = ln(Q^/Qg). Note that T depends on p = ln(l/a;) through x- 
dependence of Qs- For illustrative purposes let’s assume, following Ref. [2^ . 
that 

TipQ- p) = Qip- Pq) (7) 

representing a step-like front moving with rapdity y = ln(l/x), which plays 
a role of an evolution time. Averaging ([^ with probability distribution ([^ 
gives: 


{X(pq)) 


oo 

= jdpP(p)T(pQ- p)) 

— oo 



f\ogiQ^/Qlix))\ 

V y2a )- 


So the new scaling variable is 

logr = -\ogiQ‘^/Qlix)). 
a 


( 8 ) 

(9) 
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It is also straightforward to check that if T is a function oi pq — p then, after 
integrating over fluctuations, T —)■ (T(pp — pq)) so that geometrical scaling 
is preserved so long as the width of fluctuations is rapidity independent. 

We do not have any explicit model leading to the fluctuations with the 
constant width, although - given the phenomenological arguments above - 
such fluctuations are compatible with the DIS data. The fluctuations of the 
saturation momentum we are employing here are caused by an evolution 
with energy of the initial conditions of the forward scattering amplitude and 
are associated with the low density region (he. with a region where the 
amplitude is small, which means rather large Bjorken x’s). But for large 
x's parton distributions are non-perturbative, and the fluctuations can be of 
order one. In contrast, it was established that pomeron loops effects I2SI. 
once at small x, are small, and therefore they cause only a small change in 
the scale of rapidity fluctuations as one goes to higher and higher energies. 
Large nonperturbative fluctuations in turn, propagate to small x’s i.e. to the 
region we are interested in. 

In the following shall ignore the rapidity dependence and treat cr as a 
constant at hxed LHC energy. Obviously one needs to include beam energy 
dependence if this result would be used over a very wide range of energy from 
RHIC to the LHC. 

In the analysis of pH collisions which we present below, we will con¬ 
sider fluctuations only of the nucleon. Including the fluctuations of the 
saturation scale in the nucleus might be done, but we expect a small ef¬ 
fect. This is because the fractional change in the saturation momentum 
squared of the nucleus should be reduced by random walks. We might expect 
^Ql/Ql 1 /V^HHpart where iV^part is the number of participating nucleons 
in the target nucleus. In pp collisions, we would need to consider fluctuations 
in both protons. Without such inclusion of the effects of fluctuations of the 
nuclear saturation momentum, the results we present here should be taken 
to demonstrate qualitatively and semi-quantitatively the effects of fluctuat¬ 
ing saturation momenta, and quantitative conclusions should await a more 
systematic analysis. 

2. Limitations of the Wounded Nucleon and Color Glass Descrip¬ 
tions of the LHC pPb Data 

In the wounded nucleon (WN) model [3l H] multiplicity is given in terms of 
the number of wounded nucleons multiplied by a function describing radiation 
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of particles by one struck participant: 

dN 

— = WpFp{r]) + waFa{v)- ( 10 ) 

Total number of participants is therefore A^part = Wp + wa = 1 + wa- The 
only difference between functions Fp and Fa is that they radiate particles 
into different hemispheres. Therefore 


Fpiv) = Fiv): FAiv) = F{-r]). 


( 11 ) 


In order to test WN model against data it is convenient to introduce sym¬ 
metrized and antisymmetrized spectra |1]: 


dN 

dr] 

dN 

dr] 


sym 


asym 


= ^(>)) + ^(->)) = (», + yJA) (F{r,) + f (-,,)), (12) 

dN dN 

= ^(h) - ^(-h) = {Wp - wa) {F{v) - F{-r])). (13) 


Equation (12) implies that ratios of symmetrized spectra at different centrali¬ 


ties are independent of t] and are equal to the ratio of number of participants. 
Ratios of symmetrized spectra to the centrality class with lowest number of 
participants are plotted in the left panel of Fig. For clarity only in three 
cases error bands due to the uncertainties of iVpart are also displayed. We can 
see that an overall magnitude of these ratios agrees with the data, however 
- especially for large multiplicities - they are not r] independent in contrast 
to the prediction of the WN model. 

On the other hand antisymmetrized spectra are proportional to iVpart — 1, 
assuming that on the proton side there was only one participant. We have 
£t the slope of the antisymmetrized spectra (which are to a good accuracy 
straight lines) to the highest iVpart data, obtaining predictions for other cen¬ 
tralities. This is plotted in the right panel of Fig. for three centrality 
classes. Shaded bands around the theoretical lines correspond again to the 
uncertainties of iVpart- We see that here the wounded nucleon model fails to 
describe the data and that slopes of the antisymmetrized spectra have some 
additional dependence on iVpart- 

In the case of CGC, multiplicity is given in terms of the saturation 
scales [26]: 

dN „ . QY 


Y^ = S±Ql{2 + ln^ 


dy 


Ql 


(14) 
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Figure 3: Wounded nucleon model predictions for symmetrized (left) and antisymmetrized 
(right) spectra defined in Eqs. (|l^ and (131. In the left panel ratios of symmetrized 


distributions to the 60-90 % (c8) centrality class are plotted together with the straight 
lines corresponding to the ratios of respective A^part- Error bands due to the uncertainty 
of -/Vpart are for clarity plotted only for three centrality classes: 0-1 % (cl - top, blue), 
10-20 % (c4 - middle, orange) and 40-50 % (c7 - bottom, light blue). In the right panel 
antisymmetrized spectra (13) are plotted for three centrality classes 0-1 % (blue), 10- 
20 % (orange) and 60-90 % (yellow). Shaded bands correspond to the uncertainty of 
A^part- Experimental points are from ATLAS Collaboration. 


and the multiplicity in pA collisions is totally driven by the proton saturation 
scale. Here S± is a transverse area corresponding to an overlap of colliding 
hadrons. In the case of heavy ion and pH collisions it has rather well de- 
hned geometrical meaning, whereas in pp (small system) scattering it is a 
parameter related to the multiplicity of a given event. As such it may also 
fluctuate, however in the present paper such fluctuations are neglected. For 
hxed impact parameter 

Pa = In ^ ~ In A^part (15) 


has only logarithmic dependence on A^part- 

Formulae (14) and (15) predict both energy (denoted hereafter as IF = 
y/s) dependence and rapidity (denoted below as y) dependence and also A^part 



























dependence of multiplicities through the dependence of the saturation scales 
on these quantities [26] : 


= Ql exp(A 2 /), 

exp(-A|/) (16) 


where we take for A = 0.32 in agreement with recent analysis [27] of the 
combined DIS data from HERA. Here Wq is an arbitrary energy scale taking 
care of the proper dimension of the saturation scales. 

Pseudo-rapidity rj, which is measured, is related to rapidity y in the fol¬ 
lowing way 


2/(h) 



-|- sinh^ T] + sinh rj 
1 ]^ + sinh^ T] — sinh t] 


(17) 


The corresponding Jacobian reads: 


h{ri) = 


coshr; 


tJq + sinh^ T] 


where 


ho = 


+ Pt 
Pt 


(18) 


(19) 


In the following we shall use pg = 1.35. 

In Fig. 1^ we plot theoretical predictions of Eq. (14) and the corresponding 
data from ATLAS and ALICE. One can see that CGC formula fails to de¬ 
scribe overall normalization and the slope (note that theoretical curves have 
been normalized at p = 0). In the remainder of this paper we shall show 
that inclusion of the fluctuations of the proton saturation scale is going to 
improve upon both issues 


3. Fluctuations of the saturation scale 


In this Section we shall derive formulae for the multiplicity of gluons 
in the fluctuating case with probability distribution given by Eq. (|^. For 
fluctuating saturation scale < Q\ multiplicity takes the form given by 
Eq. (14) with Q'i replaced by and for Q\ < we simply interchange 
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Qa ^ Q- With this, we have ignored a possible slow dependence npon the 
satnration momentnm of the proton that affects the transverse area. This 
is slowly varying, like a logarithm sqnared, according to the resnlt in the 
introdnction, bnt we will later need to inclnde snch an effect to properly 
normalize the psendorapidity distribntion that we compnte. 

Now we can compnte the mnltiplicity as 


dy 


S±Ql 


CpA 


dp {pA- P + 2)epp{p) + dp (p- pa + 2)e^-^P(p) 

J PA 


( 20 ) 


The hrst integral reads 


h = 


"PA 


ip {PA- P + 2)e'’P(p) 


( 21 ) 


2/2 r 




I 

and the remaining integral is 




a 


+ 2)(1 + Erf 


Pa 


a 


72, 


a 


h = 


dp {p- Pa + 2)eP^P{p) 


( 22 ) 


'PA 

pPA 




I 


< ae 


p\Po 


+ 


-(2-p7(l-Erf 


PA 


72, 


cr 


This gives 


f = S.QJ(A + /.) 


(23) 


The satnration scales entering Eqs. (|21[)-(|23|) are given by Eq. (16); they 
depend on energy and rapidity, and - in the case of the nnclens - also on Npart- 
Note that in the limit where a —)■ 0, i.e. in the limit with no flnctnations. 


/2 —)■ 0 and Ji —)• (2 + in agreement with Eq. (14). 

In order to illnstrate the effect of flnctnations we plot in Fig. the pA 
mnltiplicity for two different centrality classes for the non-flnctnating case 
(lower bine cnrve) and for the case with flnctnations inclnded for cr = 1 
(middle orange line) and a = 2 (npper green line). We see clearly two effects 
whose strength depends on centrality: hrst, an overall normalization is in¬ 
creasing with increasing a and second, distribntions are getting hatter when 
a is rising. For small nnmber of participants we observe even a connterintn- 
itive behavior that mnltiplicity is smaller on the nnclens side (npper green 
line in the left panel of Fig. |^. 
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Figure 4: Multiplicity dependence on rapidity y for the non-fluctating case (blue) and for 
the case with fluctuations with <7 = 1 (orange) and 2 (green). Left plot corresponds to 
-^part = 6.59 and the right one to fVpart = 24.13. 


4. Multiplicity distributions in pA 

Recently ATLAS and ALICE collaborations at the LHC have published 
multiplicity distributions in pPb collisions at 5.02 TeV for different centrality 
classes summarized in Tables and Throughout this paper for ATLAS we 
have used A^part determination by means of the fluctuating Glauber-Gribov 
model with G = 0.55 and for ALIGE we use centrality classes determined 
by so called VOA method (for details see Refs. [ZlIH])- This data is shown in 
Figs. [^and|^ In Fig. |^we also plot theoretical predictions for multiplicity 
distributions without fluctuations corresponding to Eq. (14). We see that 


GGG predictions are too steep. Moreover, they are not properly normalized, 
so we have adjusted the overall factor S± by normalizing the distributions at 
r] = 0. In this way S± does - as remarked in the Introduction - depend on 
A^part, and this dependence is plotted in Fig. 

Next in Fig. we superimpose over the data theoretical predictions of 
Eqs. (21)-(23) with fluctuations included. Theoretical predictions are again 
normalized at rj = 0 and S± dependence on A^part is plotted in Fig. 

Before we discuss the effects of the fluctuations let us stress that we have 
tried to avoid playing with free parameters. Therefore we have kept A hxed to 
the DIS value of 0.32, we have kept hxed to 1.35. The only free parameter 
was a and normalization. We have found that the best description of both 
ATLAS and ALIGE data is for a = 1.55. 
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class 

% 

(-^part) 

cl 

0-1 

24.13 

c2 

1-5 

19.42 

c3 

5-10 

16.45 

c4 

10-20 

13.70 

c5 

20-30 

11.16 

c6 

30-40 

9.20 

c7 

40-60 

6.59 

c8 

60-90 

3.55 


Table 1: Mean number of participants in different centrality classes as determined by 
ATLAS by Glauber-Gribov model with fl = 0.55 


class 

% 

(-^part) 

cl 

0-5 

15.7 

c2 

5-10 

14.0 

c3 

10-20 

12.7 

c4 

20-40 

10.4 

c5 

40-60 

7.42 

c6 

60-80 

4.81 

c7 

80-100 

2.94 


Table 2: Mean number of participants in different centrality classes as determined by 
ALICE by so called VOA method. 
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Figure 5: Multiplicity dependence on pseudorapidity 77 given by formula (14) for ATLAS 
(left) and ALICE (right) for different centrality classes described in the text. Theoretical 
curves have been normalized at 77 = 0 (see Fig. [^. 


We see that fluctuations make rj distributions flatter and the A^part de¬ 
pendence of S± is also weaker than in the case with no fluctuations. The 
dependence of S± can be well approximated by the following formulae: 

SALICE ^ ^Q_Qg ^ Q_g4 in(iVp,rt))' (24) 

that take into account cross-section fluctuations with iVpart discussed in the 
Introduction. The strength of these fluctuations is, however, very different 
in the case of ATLAS and ALICE. 


5. Conclusions 

The original idea of studying multiplicity fluctuations in pA collisions 
was that they would allow probes of varying density regions inside the nu¬ 
cleus. The results we present in this paper suggest that the gluon content 
of the proton also itself changes when the multiplicity of final state particles 
varies. The analysis we present here attempts a semi-quantitative compari¬ 
son to experimental data. This data is not simply explained by the simple 
applications of the theory of the Color Glass Condensate or of the Wounded 
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Figure 6 : Multiplicity dependence on pseudo-rapidity 77 for the fluctuating case with a = 
1.55. Left plot corresponds to ATLAS whereas the right one to ALICE. Different curves 
correspond to the centrality classes defined in Tables [T] and 


Nucleon model. In both cases the shape of charge particle rapidity distribu¬ 
tion is determined by the participant number A^part, and in both cases the 
simplest models of such a dependence fail to describe the data. 

The idea that cross sections might fluctuate as an explanation for var¬ 
ious phenomena in heavy ion collisions is not new. This was used in the 
Glauber-Gribov implementation to compute the number of participants in 
the ATLAS results described above. Motivation was found for the use of hy¬ 
drodynamics in high multiplicity pA events at the LHG [28] . There have also 
been suggestions that this might be required from attempts to describe the 
ridge in pA collisions [22] • Furthermore fluctuations have been implemented 
phenomenologically in a Monte Garlo study of pA collisions done in Ref. [32] 
where the authors consider multiplicity variations for hxed Apart given by the 
negative binomial distribution (NBD) that is obtained in the glasma model 
[311 [ 32 ] . There our formula (14) corresponds to the average number of parti¬ 
cles entering NBD. It is, however, difficult to assess the energy dependence of 
these fluctuations as it stems from the energy dependence of the NBD itself 
and of the impact parameter prohle that selects given Apart (that in a sense 
would correspond to the fluctuations of the nucleus saturation scale). 

On the other hand, the arguments for a fluctuating saturation momen¬ 
tum have been primarily theoretical [IS]^[2S]- Interest in such fluctuations 
diminished somewhat when no evidence of diffusive scaling has been found 
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Figure 7: Normalization of multiplicity distributions shown in Fig[^as functions of A^part- 
Left panel is for ATLAS and right panel for ALICE. Upper (orange) points correspond 
to the non-fluctuating case, whereas lower (blue) points correspond to fluctuations with 
(T = 1.55. Parametrization of Eq. (24) is shown as a solid (blue) line. 


in DIS [201 122] and after it was argued that the evolution of the strength 
of such fluctuations was very weak with varying rapidity. However, if the 
fluctuations are largely rapidity independent - as it is the case of the present 
paper - geometrical scaling variable r is not modihed and there is no con¬ 
tradiction between our results and the absence of diffusive scaling in DIS. 
Another point is the energy dependence of such fluctuations. This needs to 
be tested in comparative studies at both RHIC and LHC energies. 

Our hnding that the fluctuations can account for the data must be taken 
with much caution. It is simply a hrst attempt and the number extracted 
for the width of the fluctuations of the logarithm is rather large, and exceeds 
by a factor of 1.5 the one of recent analysis of Ref. [23] • This may be be¬ 
cause in our simple analysis, we do not consider fluctuations of the saturation 
momentum of the nucleus. It is not obvious how to do this, since the fluc¬ 
tuations would be local over some transverse scale size of order the inverse 
saturation momentum. Similarly, we expect there to be local structure in the 
proton since the size scale of the saturation momentum is small compared to 
the proton size. The issue of a local fluctuating saturation momentum was 
addressed by lancu and McLerran [33] within the context of conformal held 
theory, but there is as yet no practical implementation of such ideas useful 
for understanding experimental data. Furthermore, in our simple analysis 
we have neglected other sources of huctuations that have been taken into ac¬ 
count in Ref. [23] • Finally, since we do not include any dumping for large y's 
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(as has been done e.g. in Ref. [2S] where large x suppression for gluon den¬ 
sities has been included), our calculation is most reliable in the mid rapidity 
region. Evidently a more realistic model must contain such effects. 

If we accept the idea of a fluctuating saturation momenta as a working 
hypothesis, there are a variety of tests one can imagine. Fluctuations in the 
multiplicity for very high multiplicity predicted in the IP-glasma model will 
be modified by such effects [35l |33] . The correlation between multiplicity and 
jet production will also be changed, since higher saturation momenta in the 
proton means more gluons and therefore greater probability of producing a 
jet. 

The original idea of studying multiplicity fluctuations in pA collisions was 
that they would allow probes of higher density regions inside the nucleus. The 
results we present in this paper suggest that the gluon content of the proton 
also itself changes when the multiplicity of hnal state particles varies. 
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